Kinetic approach to the cluster liquid-gas transition 
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The liquid-gas transition in free atomic clusters is investigated theoretically based on simple 
unimolecular rate theories and assuming sequential evaporations. A kinetic Monte Carlo scheme is 
used to compute the time-dependent properties of clusters undergoing multiple dissociations, and 
two possible definitions of the boiling point are proposed, relying on the cluster or gas temperature. 
This numerical approach is supported by molecular dynamics simulations of clusters made of sodium 
atoms or Ceo molecules, as well as simplified rate equation. 
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Boiling is often considered as archetypal of first-order 
phase transitions in bulk matter. Recent experimental 
evidence 0, E3i suggests that the liquid-gas tran- 
sition also occurs in finite atomic systems. The caloric 
curves measured for sodium Q, hydrogen Q or stron- 
tium Q clusters exhibit a plateau or a backbending, be- 
lieved to be signatures of a phase transition rounded by 
size effects. Similar conclusions have been inferred from 
collisions between gold nuclei The experiments per- 
formed by Pochodzalla and coworkers have since moti- 
vated a significant amount of theoretical work to help the 
search for apossible equation of state for nuclear matter 

For more than two decades, fragmentation has been 
recognized as one of the most convenient ways of ac- 
cessing fundamental cluster properties such as binding 
energies or temperatures. These quantities were related 
to each other through calorimetric experiments on the 
solid- liquid phase change [111 Il2j | . The crucial role of the 
observation time was also noticed at an early stage, and 
more fully understood by Klots 0] who introduced the 
concept of the evaporative ensemble. Free clusters are 
never strictly stable when their energy exceeds a certain 
dissociation threshold, as fragmentation will occur, pos- 
sibly very late. Measurements on free clusters are thus 
conducted on species resulting from hotter and bigger 
clusters, and which have sufficiently cooled down so as 
not to evaporate further. Unfortunately, the relatively 
long times involved in experiments (usually the /is-ms 
range) have prevented Klots' ideas from being exploited 
in later theoretical studies. While melting can be con- 
veniently addressed using Monte Carlo (MC) or molecu- 
lar dynamics (MD) simulations by artificially keeping the 
cluster in a container [LH , fragmentation of a cluster into 
vacuum is an out-of-equilibrium phenomenon. The dif- 
ficulty of accounting for the time variable has generally 
been circumvented by considering lattice-gas or percola- 
tion models 0, as well as periodic boundary condi- 
tions 

Unimolecular rate theories provide a general frame- 
work to describe single dissociation events accurately and 
over long time scales Here we introduce a kinetic 



Monte Carlo (KMC) scheme based on such rate theories 
to calculate the caloric curves of clusters across the liquid- 
gas transition for arbitrarily long observation times. 

We start by illustrating the boiling problem by show- 
ing the results of classical MD simulations carried out 
on three selected clusters, namely Na55, (C6o)4o> an d 
a model binary cluster X13Y42. These systems are de- 
scribed using an explicit many-body potential [lfij . the 
Girifalco potential ^21 > an d Lennard- Jones (LJ) interac- 
tions, respectively. For the LJ cluster, we chose £xx = 1, 
£yy = 1/2 and £xy = l/v 7 ^, all distance units and 
masses being set to one. For each cluster, MD trajectories 
are performed starting from the lowest-energy minimum, 
at increasing total energies and zero angular and linear 
momenta. After some observation time t we determine 
the biggest remaining cluster as the largest set of con- 
nected atoms or molecules, two atoms being connected 
when their distance is smaller than a cut-off value r cu t 
chosen as twice the equilibrium distance. The instanta- 
neous temperature of this cluster is then calculated after 
removing the overall translation and rotation contribu- 
tions. The average cluster size and temperature obtained 
over 100 independent trajectories (20 for the longest wait- 
ing times) are represented for the three clusters in Fig. ^ 

All caloric curves exhibit two clearly distinct regimes. 
At low energies, the roughly linear increase of tempera- 
ture characterizes the condensed phases as solidlike and 
liquidlike. For the sodium and, to a lesser extent, for the 
LJ cluster, a small inflection marks the onset of melt- 
ing. At high energies, the cluster temperature drops and 
reaches a plateau. The plateau temperature is lower and 
the change between the two regimes is sharper as the 
waiting time increases. This behavior is consistent with 
the observed smaller average cluster size, and results from 
a stronger evaporative cooling. In the cluster of fullerene 
molecules, no evidence for the melting transition is seen 
on the caloric curves, in agreement with the lack of sta- 
bility of the liquid phase known in bulk Cgo 18] . The 
backbending seen for the shortest waiting time (0.5 ns) 
was confirmed by another set of 100 independent MD 
trajectories. However, without any container preventing 
dissociation, it is hard to attribute this feature to melting 



2 



1250 




Energy (eV/atom) Energy (eV/molecule) 
1250 




Energy (eV/molecule) Energy (eV/molecule) 




Energy (e) Energy (e) 



FIG. 1: Caloric curves and average cluster size after various 
waiting times in MD simulations, (a) and (b) Nass; (c) and 
(d) (Ceo)4o; (e) and (f) binary LJ cluster X13Y42. Cluster 
sizes are given relative to the initial size. 

or to boiling. The binary LJ cluster was constructed as 
a two-layer icosahedron, with the most strongly bound 
atoms X in the center. Provided that the waiting time is 
long enough, X atoms may dissociate, but only after the 
Y atoms have evaporated. The fragmentation of X atoms 
is shown on the caloric curve corresponding to t = 10 ns 
by the slight increase at high excitation energies. 

The time-dependence of fragmentation caloric curves 
can also be calculated from simple cluster models, sim- 
ilar to those introduced by Bixon and Jortner for the 
isomerization problem [l^. Our first assumption is that 
clusters are heated adiabatically, allowing fragmentation 
to occur through sequential loss of monomers HI. We 
use simple rate theories to describe each dissociation 
step. The dissociation rate k n {E n ) of cluster X„ into 
X„_i depends on the dissociation energy A„ and the 
total energy E n through the harmonic RRK approach 
HH], namely k n (E n ) = v (l - A n /E n ) 3n - 6 . This well 
known formula contains a single time scale parameter 
Uq, related to the typical vibrational period. Once evap- 
oration has occured, the cluster loses a part of its in- 
ternal energy, which is described more accurately us- 



ing the Weisskopf theory |22j. The probability distri- 
bution of the kinetic energy released (KER), p n (e,E n ), 
is p n (e,E n ) oc e(E n — A„ — e) 3n ~ 7 . The dissociation 
of rigid molecules would be described similarly, changing 
the factor e and the number 3n — 7 of degrees of freedom 
to e 3 / 2 and 5n — 7 for linear molecules, and to e 2 and 
6rt — 7 for tridimensional molecules, respectively. 

The problem can be further simplified by assuming 
that the lifetime of the parent cluster is l/k n and that 
evaporative cooling removes the average KER (e)„ = 
2(E — A„)/(3n — 7), leading to the energy of the product 
X„_i, E n -i — E n — A„ — (e) n . The multifragmenta- 
tion problem is thus reduced to computing all successive 
dissociation rates, from which the survival probabilities 
p n (t) of all cluster sizes are calculated following a master 
rate equation dp n /dt = k n p n — fc n _ip n _i, solved exactly. 
A similar "mean-field" technique was used by Hervieux et 
al. to calculate branching ratios in the collision-induced 
fragmentation of Na^J" [23j . 

Beyond this approximate treatment, the kinetic Monte 
Carlo method |24j accounts for the continuous character 
of dissociation times and KER distributions. Starting 
with the parent cluster size n at total energy E n , the 
evaporation rate k n (E n ) is calculated and the KER e 
is chosen randomly from the distribution p n (e, E n ). The 
energy is decreased by A„ + e, and the time is updated by 
the quantity — log(ran)//c„, where < ran < 1 is a ran- 
dom real number. This process is iterated until cither the 
waiting time has been exceeded or the remaining energy 
is below the next dissociation limit. The final cluster has 
k atoms and the energy Ek, its temperature is obtained 
in the harmonic approximation ksTk = Ek/(3k — 6). 

The KMC and mean-field methods have been com- 
pared on the simplest case of the multiple dissociation 
of a 100-particle model cluster with A„ = 1 for all n. 
Figs. Eta) and (b) show the caloric curve and the av- 
erage final cluster size obtained from the KMC calcula- 
tions. The mean field values were not reported, as they 
are undistinguishable from the stochastic data. It should 
be noted that these results are not significantly affected 
when the energetics of dissociation are described using 
RRK theory, p n (e, E n ) cx {E n - A„ - e) 3 "" 6 . 

The MD results of Fig. ^ are well reproduced by the 
predictions of this simple model, especially the plateau 
of the vapor phase and the increasing sharpness of the 
transition for long observation times. The absence of any 
feature on the caloric curve for t — > demonstrates the 
kinetic character of the cluster liquid-gas transition. On 
the other hand the asymptotic curve for t — > 00 exhibits 
multiple tiny backbendings corresponding to consecutive 
dissociations leaving the cluster perfectly cold. We have 
simulated the effect of a more stable cluster at size n* = 
90 by imposing A„» = 1.5. The influence of this magic 
cluster on the curves of Figs.EIc) and (d) is rather local. 
The greater stability of the X„. cluster induces a delay 
on the subsequent dissociations, but the overall plateau 
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FIG. 2: Caloric curves and average final cluster size obtained 
from KMC calculations for model clusters without (a and b) 
or with (c and d) a specially stable fragment at n = 90. The 
waiting times are given in units of 1/Vo, and the cluster sizes 
are relative to the initial size 100. 



temperatures do not vary much. 

The kinetic Monte Carlo method can be used in more 
general situations, for which the rate equation approach 
is not practical. Clusters with competing fragmentation 
channels turn into products having different internal en- 
ergies depending on their history. For instance, hetero- 
geneous clusters X„Y p can dissociate in multiple ways: 

/ \ 

/ \ S \ 

X„_ 1 Y p _ 1 

Statistically the X„_iY p _i fragment does not have the 
same energy depending on whether it was produced from 
X n _!Y p or from X„Y p _!. For the above example, and 
within the Weisskopf theory do the two average energies 
of the X„_iY p _i product differ by 



A£„-i, P -i = 2(A X -A Y ) 

3(n+p) - 12 



[3(n + p) - 7] • [3(n+j>) - 10] 



, (1) 



depending on which of X or Y was emitted first. Only 
at large sizes or for similar binding energies Ax — Ay 
the two fragmentation routes become equivalent. More- 
over, clusters may be prepared at fixed temperature 
rather than fixed energy. The initial energy is then ran- 
domly picked from the canonical distribution p n (E) oc 
E 3n - 6 exp(-E/k B T) as the first Monte Carlo step. 



By mimicking the actual dissociation cascade followed 
by each individual cluster, the KMC technique offers a re- 
alistic statistical description of multifragmentation. The 
caloric curves of model binary clusters X„Yioo-n were 
calculated assuming that X and Y atoms are bound to 
the clusters by the energies Ax = 1 and Ay = 1/2, re- 
spectively. The caloric curves represented in Fig.^a) for 
several compositions exhibit two steps, each associated 
with the boiling of Y, then X particles. This two-step 
liquid-gas transition is somewhat similar to the multi- 
step melting process often seen in cluster simulations |25j . 
As more Y atoms are added, a larger energy is required 
to evaporate them, therefore the gas transition of the 
remaining X cluster is initiated at higher temperatures. 
The different onsets of boiling for X and Y particles in 
Fig- Wi c ) an d (d) confirm this two-step transition. 
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FIG. 3: Final cluster (a) and gas (b) temperatures versus ini- 
tial temperature in the boiling transition of model binary clus- 
ters with various compositions, for the waiting time 10 6 /Vo- 

We now address possible definitions of the boiling tem- 
perature, based on the previous results. The onset of the 
drop in the slope of the caloric curve provides one possible 
estimate. However, especially at short waiting times, the 
variations of the average remaining cluster size seem to 
be more reliable. Another observable consists of looking 
at the kinetic energy of the evaporated atoms, in a fash- 
ion somewhat more closely related to experimental con- 
ditions. Using the same notation as above for the initial 
and final cluster sizes and energies, and assuming that all 
n — k evaporated atoms behave like a perfect gas, a gas 
temperature can be defined by 3(ji-fc)fcgT/2 = E n —E^. 

The variations of the average gas temperature across 
the boiling transition are shown in Fig-Eld) for the model 
binary clusters. Boiling is made evident by the sudden 
increase of the gas temperature, which provides us with 
another characterization of the liquid-gas point. The sec- 
ond boiling point involving the more strongly bound X 
atoms does not lead to significant variations of the gas 
temperature, because X and Y atoms contribute equally 
once in the gas phase. A better estimate of the liquid-gas 
transition temperature of the cluster of X atoms is found 
with the variations of its final size [see Fig. IHd)] . 

The three aforementioned definitions of the boiling 
point yield similar values in a broad range of situations 
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FIG. 4: Boiling temperature of model 100-atom homogeneous 
cluster versus waiting time. Inset: boiling temperatures of the 
X and Y parts of model X n Yioo-n clusters versus composition 
n, for the waiting time t = 10 6 /fo- 



covered by our simple models. The dependence of the 
observation time on the boiling temperature of the ho- 
mogeneous cluster is represented in Fig. 0] As expected, 
waiting longer favors evaporation, hence a lower boiling 
temperature. In most experiments 0,0 0,13] the waiting 
time is at least 5 orders of magnitude larger than the typ- 
ical vibrational period. Fig. 0] suggests that the results of 
these experiments should remain stable by less than 5% 
if the waiting time is doubled or halved. However, similar 
measurements on trapped clusters might show some de- 
viation since the trapping time may exceed seconds |26| . 

In Fig. 0| we also show how the presence of more weakly 
bound (Y) atoms influence the boiling point of an (X) 
cluster. The fragmentation temperature of Y atoms does 
not significantly change with their initial composition. 
However, as more Y atoms are added, the X product gets 
progressively colder due to the more numerous evapora- 
tions. This explains why the boiling temperature of the 
X cluster decreases with its composition. 

The KMC technique outlined in this paper could be 
improved using more accurate unimolecular rate theo- 
ries, such as Phase Space Theory, which incorporates an- 
harmonic densities of states as well as a rig orous treat- 
ment of angular momentum constraints The main 
assumption of the present approach is that boiling in 
slowly heated clusters occurs sequentially. Our simu- 
lation results seem to validate this hypothesis. As in 
glasses (2?J , the arbitrarily long time scales reached by 
the present statistical approach make it a useful alterna- 
tive to molecular dynamics. We believe that it forms a 
bridge between unimolecular rate descriptions and multi- 
fragmentation models that assume thermal equilibrium. 

As a first application of the KMC method to the clus- 
ter dissociation problem, we have considered the liquid- 
gas transition. Our results emphasize the important role 
played by the observation time on the caloric curves. 



They also indicate that multiple-step boiling transitions 
could be detected in heterogeneous clusters. Beyond 
these examples, a more complete interpretation of recent 
experiments 0, 0, especially on clusters exhibiting 
competing dissociation channels could be anticipated. 
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